Coop Data

NOAA/NWS Cooperative Observer Network

Downloading from the Iowa Environmental Mesonet. (Data taken from the neighboring COOP analysis

Station: DEL NORTE 2E, CO2184

COOP_stations <- read.csv("C:/Users/13074/Documents/ESS580/thesis_project/All_COOP_stations/data_raw/COOP_stations.csv", header = TRUE)
#str(COOP_stations)
COOP_stations$Date <- ymd(COOP_stations$day)
del_norte_clean <- COOP_stations %>% # filter for the timeframe
  filter(station == "CO2184") %>% 
  addWaterYear() %>%
  mutate(daymonth = format(as.Date(Date), "%d-%m")) %>% 
  group_by(waterYear)%>% 
  mutate(waterDay = (as.integer(difftime(Date, ymd(paste0(waterYear - 1 ,'-09-30')), units = "days")))) %>%   na.omit()
del_norte_clean <- del_norte_clean %>% 
  mutate(avg_T_c = (highc+lowc)/2) %>% 
  filter(waterYear <= 2022)
write.csv(del_norte_clean,"C:/Users/13074/Documents/ESS580/thesis_project/del_norte/data_clean/del_norte_clean.csv", row.names = FALSE)

Figure check

str(del_norte_clean)
## grouped_df [47,389 x 13] (S3: grouped_df/tbl_df/tbl/data.frame)
##  $ station     : chr [1:47389] "CO2184" "CO2184" "CO2184" "CO2184" ...
##  $ station_name: chr [1:47389] "DEL NORTE 2E" "DEL NORTE 2E" "DEL NORTE 2E" "DEL NORTE 2E" ...
##  $ lat         : num [1:47389] 37.7 37.7 37.7 37.7 37.7 ...
##  $ lon         : num [1:47389] -106 -106 -106 -106 -106 ...
##  $ day         : chr [1:47389] "1893/01/01" "1893/01/02" "1893/01/03" "1893/01/04" ...
##  $ doy         : int [1:47389] 1 2 3 4 5 6 7 8 9 10 ...
##  $ highc       : num [1:47389] 2.2 1.7 9.4 7.2 7.2 7.8 7.2 5 4.4 7.2 ...
##  $ lowc        : num [1:47389] -13.3 -12.2 -11.7 -11.7 -12.8 -13.3 -13.3 -12.2 -11.1 -11.7 ...
##  $ Date        : Date[1:47389], format: "1893-01-01" "1893-01-02" ...
##  $ waterYear   : num [1:47389] 1893 1893 1893 1893 1893 ...
##  $ daymonth    : chr [1:47389] "01-01" "02-01" "03-01" "04-01" ...
##  $ waterDay    : int [1:47389] 93 94 95 96 97 98 99 100 101 102 ...
##  $ avg_T_c     : num [1:47389] -5.55 -5.25 -1.15 -2.25 -2.8 -2.75 -3.05 -3.6 -3.35 -2.25 ...
##  - attr(*, "groups")= tibble [130 x 2] (S3: tbl_df/tbl/data.frame)
##   ..$ waterYear: num [1:130] 1893 1894 1895 1896 1897 ...
##   ..$ .rows    : list<int> [1:130] 
##   .. ..$ : int [1:273] 1 2 3 4 5 6 7 8 9 10 ...
##   .. ..$ : int [1:365] 274 275 276 277 278 279 280 281 282 283 ...
##   .. ..$ : int [1:365] 639 640 641 642 643 644 645 646 647 648 ...
##   .. ..$ : int [1:366] 1004 1005 1006 1007 1008 1009 1010 1011 1012 1013 ...
##   .. ..$ : int [1:365] 1370 1371 1372 1373 1374 1375 1376 1377 1378 1379 ...
##   .. ..$ : int [1:365] 1735 1736 1737 1738 1739 1740 1741 1742 1743 1744 ...
##   .. ..$ : int [1:365] 2100 2101 2102 2103 2104 2105 2106 2107 2108 2109 ...
##   .. ..$ : int [1:365] 2465 2466 2467 2468 2469 2470 2471 2472 2473 2474 ...
##   .. ..$ : int [1:365] 2830 2831 2832 2833 2834 2835 2836 2837 2838 2839 ...
##   .. ..$ : int [1:365] 3195 3196 3197 3198 3199 3200 3201 3202 3203 3204 ...
##   .. ..$ : int [1:365] 3560 3561 3562 3563 3564 3565 3566 3567 3568 3569 ...
##   .. ..$ : int [1:366] 3925 3926 3927 3928 3929 3930 3931 3932 3933 3934 ...
##   .. ..$ : int [1:365] 4291 4292 4293 4294 4295 4296 4297 4298 4299 4300 ...
##   .. ..$ : int [1:365] 4656 4657 4658 4659 4660 4661 4662 4663 4664 4665 ...
##   .. ..$ : int [1:365] 5021 5022 5023 5024 5025 5026 5027 5028 5029 5030 ...
##   .. ..$ : int [1:366] 5386 5387 5388 5389 5390 5391 5392 5393 5394 5395 ...
##   .. ..$ : int [1:365] 5752 5753 5754 5755 5756 5757 5758 5759 5760 5761 ...
##   .. ..$ : int [1:365] 6117 6118 6119 6120 6121 6122 6123 6124 6125 6126 ...
##   .. ..$ : int [1:365] 6482 6483 6484 6485 6486 6487 6488 6489 6490 6491 ...
##   .. ..$ : int [1:366] 6847 6848 6849 6850 6851 6852 6853 6854 6855 6856 ...
##   .. ..$ : int [1:365] 7213 7214 7215 7216 7217 7218 7219 7220 7221 7222 ...
##   .. ..$ : int [1:365] 7578 7579 7580 7581 7582 7583 7584 7585 7586 7587 ...
##   .. ..$ : int [1:365] 7943 7944 7945 7946 7947 7948 7949 7950 7951 7952 ...
##   .. ..$ : int [1:366] 8308 8309 8310 8311 8312 8313 8314 8315 8316 8317 ...
##   .. ..$ : int [1:365] 8674 8675 8676 8677 8678 8679 8680 8681 8682 8683 ...
##   .. ..$ : int [1:365] 9039 9040 9041 9042 9043 9044 9045 9046 9047 9048 ...
##   .. ..$ : int [1:365] 9404 9405 9406 9407 9408 9409 9410 9411 9412 9413 ...
##   .. ..$ : int [1:366] 9769 9770 9771 9772 9773 9774 9775 9776 9777 9778 ...
##   .. ..$ : int [1:365] 10135 10136 10137 10138 10139 10140 10141 10142 10143 10144 ...
##   .. ..$ : int [1:365] 10500 10501 10502 10503 10504 10505 10506 10507 10508 10509 ...
##   .. ..$ : int [1:365] 10865 10866 10867 10868 10869 10870 10871 10872 10873 10874 ...
##   .. ..$ : int [1:366] 11230 11231 11232 11233 11234 11235 11236 11237 11238 11239 ...
##   .. ..$ : int [1:365] 11596 11597 11598 11599 11600 11601 11602 11603 11604 11605 ...
##   .. ..$ : int [1:365] 11961 11962 11963 11964 11965 11966 11967 11968 11969 11970 ...
##   .. ..$ : int [1:365] 12326 12327 12328 12329 12330 12331 12332 12333 12334 12335 ...
##   .. ..$ : int [1:366] 12691 12692 12693 12694 12695 12696 12697 12698 12699 12700 ...
##   .. ..$ : int [1:365] 13057 13058 13059 13060 13061 13062 13063 13064 13065 13066 ...
##   .. ..$ : int [1:365] 13422 13423 13424 13425 13426 13427 13428 13429 13430 13431 ...
##   .. ..$ : int [1:365] 13787 13788 13789 13790 13791 13792 13793 13794 13795 13796 ...
##   .. ..$ : int [1:366] 14152 14153 14154 14155 14156 14157 14158 14159 14160 14161 ...
##   .. ..$ : int [1:365] 14518 14519 14520 14521 14522 14523 14524 14525 14526 14527 ...
##   .. ..$ : int [1:365] 14883 14884 14885 14886 14887 14888 14889 14890 14891 14892 ...
##   .. ..$ : int [1:365] 15248 15249 15250 15251 15252 15253 15254 15255 15256 15257 ...
##   .. ..$ : int [1:366] 15613 15614 15615 15616 15617 15618 15619 15620 15621 15622 ...
##   .. ..$ : int [1:365] 15979 15980 15981 15982 15983 15984 15985 15986 15987 15988 ...
##   .. ..$ : int [1:365] 16344 16345 16346 16347 16348 16349 16350 16351 16352 16353 ...
##   .. ..$ : int [1:365] 16709 16710 16711 16712 16713 16714 16715 16716 16717 16718 ...
##   .. ..$ : int [1:366] 17074 17075 17076 17077 17078 17079 17080 17081 17082 17083 ...
##   .. ..$ : int [1:365] 17440 17441 17442 17443 17444 17445 17446 17447 17448 17449 ...
##   .. ..$ : int [1:365] 17805 17806 17807 17808 17809 17810 17811 17812 17813 17814 ...
##   .. ..$ : int [1:365] 18170 18171 18172 18173 18174 18175 18176 18177 18178 18179 ...
##   .. ..$ : int [1:366] 18535 18536 18537 18538 18539 18540 18541 18542 18543 18544 ...
##   .. ..$ : int [1:365] 18901 18902 18903 18904 18905 18906 18907 18908 18909 18910 ...
##   .. ..$ : int [1:365] 19266 19267 19268 19269 19270 19271 19272 19273 19274 19275 ...
##   .. ..$ : int [1:365] 19631 19632 19633 19634 19635 19636 19637 19638 19639 19640 ...
##   .. ..$ : int [1:366] 19996 19997 19998 19999 20000 20001 20002 20003 20004 20005 ...
##   .. ..$ : int [1:365] 20362 20363 20364 20365 20366 20367 20368 20369 20370 20371 ...
##   .. ..$ : int [1:365] 20727 20728 20729 20730 20731 20732 20733 20734 20735 20736 ...
##   .. ..$ : int [1:365] 21092 21093 21094 21095 21096 21097 21098 21099 21100 21101 ...
##   .. ..$ : int [1:366] 21457 21458 21459 21460 21461 21462 21463 21464 21465 21466 ...
##   .. ..$ : int [1:365] 21823 21824 21825 21826 21827 21828 21829 21830 21831 21832 ...
##   .. ..$ : int [1:365] 22188 22189 22190 22191 22192 22193 22194 22195 22196 22197 ...
##   .. ..$ : int [1:365] 22553 22554 22555 22556 22557 22558 22559 22560 22561 22562 ...
##   .. ..$ : int [1:366] 22918 22919 22920 22921 22922 22923 22924 22925 22926 22927 ...
##   .. ..$ : int [1:365] 23284 23285 23286 23287 23288 23289 23290 23291 23292 23293 ...
##   .. ..$ : int [1:365] 23649 23650 23651 23652 23653 23654 23655 23656 23657 23658 ...
##   .. ..$ : int [1:365] 24014 24015 24016 24017 24018 24019 24020 24021 24022 24023 ...
##   .. ..$ : int [1:366] 24379 24380 24381 24382 24383 24384 24385 24386 24387 24388 ...
##   .. ..$ : int [1:365] 24745 24746 24747 24748 24749 24750 24751 24752 24753 24754 ...
##   .. ..$ : int [1:365] 25110 25111 25112 25113 25114 25115 25116 25117 25118 25119 ...
##   .. ..$ : int [1:365] 25475 25476 25477 25478 25479 25480 25481 25482 25483 25484 ...
##   .. ..$ : int [1:366] 25840 25841 25842 25843 25844 25845 25846 25847 25848 25849 ...
##   .. ..$ : int [1:365] 26206 26207 26208 26209 26210 26211 26212 26213 26214 26215 ...
##   .. ..$ : int [1:365] 26571 26572 26573 26574 26575 26576 26577 26578 26579 26580 ...
##   .. ..$ : int [1:365] 26936 26937 26938 26939 26940 26941 26942 26943 26944 26945 ...
##   .. ..$ : int [1:366] 27301 27302 27303 27304 27305 27306 27307 27308 27309 27310 ...
##   .. ..$ : int [1:365] 27667 27668 27669 27670 27671 27672 27673 27674 27675 27676 ...
##   .. ..$ : int [1:365] 28032 28033 28034 28035 28036 28037 28038 28039 28040 28041 ...
##   .. ..$ : int [1:365] 28397 28398 28399 28400 28401 28402 28403 28404 28405 28406 ...
##   .. ..$ : int [1:366] 28762 28763 28764 28765 28766 28767 28768 28769 28770 28771 ...
##   .. ..$ : int [1:365] 29128 29129 29130 29131 29132 29133 29134 29135 29136 29137 ...
##   .. ..$ : int [1:365] 29493 29494 29495 29496 29497 29498 29499 29500 29501 29502 ...
##   .. ..$ : int [1:365] 29858 29859 29860 29861 29862 29863 29864 29865 29866 29867 ...
##   .. ..$ : int [1:366] 30223 30224 30225 30226 30227 30228 30229 30230 30231 30232 ...
##   .. ..$ : int [1:365] 30589 30590 30591 30592 30593 30594 30595 30596 30597 30598 ...
##   .. ..$ : int [1:365] 30954 30955 30956 30957 30958 30959 30960 30961 30962 30963 ...
##   .. ..$ : int [1:365] 31319 31320 31321 31322 31323 31324 31325 31326 31327 31328 ...
##   .. ..$ : int [1:366] 31684 31685 31686 31687 31688 31689 31690 31691 31692 31693 ...
##   .. ..$ : int [1:365] 32050 32051 32052 32053 32054 32055 32056 32057 32058 32059 ...
##   .. ..$ : int [1:365] 32415 32416 32417 32418 32419 32420 32421 32422 32423 32424 ...
##   .. ..$ : int [1:365] 32780 32781 32782 32783 32784 32785 32786 32787 32788 32789 ...
##   .. ..$ : int [1:366] 33145 33146 33147 33148 33149 33150 33151 33152 33153 33154 ...
##   .. ..$ : int [1:365] 33511 33512 33513 33514 33515 33516 33517 33518 33519 33520 ...
##   .. ..$ : int [1:365] 33876 33877 33878 33879 33880 33881 33882 33883 33884 33885 ...
##   .. ..$ : int [1:365] 34241 34242 34243 34244 34245 34246 34247 34248 34249 34250 ...
##   .. ..$ : int [1:366] 34606 34607 34608 34609 34610 34611 34612 34613 34614 34615 ...
##   .. ..$ : int [1:365] 34972 34973 34974 34975 34976 34977 34978 34979 34980 34981 ...
##   .. ..$ : int [1:365] 35337 35338 35339 35340 35341 35342 35343 35344 35345 35346 ...
##   .. ..$ : int [1:365] 35702 35703 35704 35705 35706 35707 35708 35709 35710 35711 ...
##   .. .. [list output truncated]
##   .. ..@ ptype: int(0) 
##   ..- attr(*, ".drop")= logi TRUE
ggplot(del_norte_clean, aes(x = Date, y = avg_T_c)) +
  geom_line() + #lwd = 2) +
  theme_few() +
  geom_smooth(method = "lm", se=FALSE) +
  ylab('Daily temperature (°C)') + 
  xlab('Date')

#Check for outliers....

#dygraph

del_norte_temp_xts <- xts(del_norte_clean$avg_T_c, order.by = del_norte_clean$Date)

dygraph(del_norte_temp_xts) %>%
  dyAxis("y", label = "Daily temperature (°C)")

Detrending Data

There are two pieces of detrended: shifting by the mean and the sigmoidal detrending.

#SF figured out the yearly average by water year

#average water year temperature

del_norte_yearly_wy_aver <- del_norte_clean %>% 
  group_by(waterYear) %>% 
  mutate(aver_ann_temp = mean(avg_T_c))
#Average temperature by day for all water years:

del_norte_daily_wy_aver <- del_norte_yearly_wy_aver %>% 
  group_by(daymonth) %>% 
  mutate(aver_day_temp = mean(avg_T_c))

#average mean temperature by day for the period of record:

del_norte_daily_wy_aver <- del_norte_daily_wy_aver %>% 
  group_by(daymonth) %>% 
  mutate(all_ave_temp = mean(del_norte_daily_wy_aver$aver_day_temp))

#str(daily_wy_aver)
# try to show all years as means. 
del_norte_daily_wy_aver2 <- del_norte_daily_wy_aver %>% 
  #filter(waterYear == "1987" | waterYear == "2021") %>%
  group_by(waterDay) %>%
  mutate(date_temp = mean(avg_T_c))
  
del_norte_daily_wy_aver2$date_temp <- signif(del_norte_daily_wy_aver2$date_temp,3) #reduce the sig figs

ggplot(del_norte_daily_wy_aver2, aes(x = waterDay, y = date_temp))+#, color = waterYear)) +
  geom_line(size= 0.7) +
  #geom_line(aes) +
  theme_few() +
  #geom_smooth(method = "lm", se=FALSE) +
  ylab('Average Daily temperature (°C)') + 
  xlab('Day of water year')

Standard Deviation

To figure out the standard deviation for each year, I want the “residual” for each daily value.

The standard deviation will be the daily residual minus the mean of the residuals by water year, summed and squared, then divided by the number of observations minus one. The square root of the resulting value of which is thus the standard deviation for the water year.

Determining residuals

del_norte_standard_dev <- del_norte_daily_wy_aver %>% 
  group_by(waterYear) %>% 
  mutate(residual = (all_ave_temp-aver_ann_temp)+avg_T_c-aver_day_temp) %>% 
  mutate(deviation = abs(residual-lag(residual)))

mean(del_norte_standard_dev$residual)
## [1] -6.725794e-17

The mean of the residuals is close enough to zero

Calculating standard deviation for the timeseries

del_norte_standard_dev_all <- del_norte_standard_dev %>% 
  group_by(waterYear) %>% 
  mutate(nmbr = n())

del_norte_standard_dev_all <- del_norte_standard_dev_all %>% 
  group_by(waterYear) %>% 
  mutate(resid_mean = mean(residual)) %>%
  mutate(sd_1 = residual-resid_mean) %>% 
  mutate(sd_2 = (((sum((sd_1)^2))/((nmbr-1))))^(0.5)) %>%
  distinct(sd_2, .keep_all = TRUE) %>% 
   select(waterYear, sd_2)

del_norte_standard_dev_all %>% 
  kable(.,'html') %>%
  kable_styling() %>%
  scroll_box(width='250px',height='500px')
waterYear sd_2
1893 2.740958
1894 3.037105
1895 3.046082
1896 2.937365
1897 3.070332
1898 3.111323
1899 3.164378
1900 2.799577
1901 3.151202
1902 3.328052
1903 3.214077
1904 2.948437
1905 2.751700
1906 2.665888
1907 3.615727
1908 2.450495
1909 3.278400
1910 3.713749
1911 3.208232
1912 2.855865
1913 2.945180
1914 2.921295
1915 2.524662
1916 3.028834
1917 3.104708
1918 3.221143
1919 3.487615
1920 3.284928
1921 3.217066
1922 3.563001
1923 2.652857
1924 3.245609
1925 3.417730
1926 2.866595
1927 3.589082
1928 3.617920
1929 2.782898
1930 3.668425
1931 2.869577
1932 3.801319
1933 3.736147
1934 2.779587
1935 3.457251
1936 2.943620
1937 3.393346
1938 3.051067
1939 3.593912
1940 3.117468
1941 3.164931
1942 3.151584
1943 3.161053
1944 3.157730
1945 2.917152
1946 3.149310
1947 3.393296
1948 3.916021
1949 3.328426
1950 3.412837
1951 3.458726
1952 3.375300
1953 4.057674
1954 3.392285
1955 3.340121
1956 3.633411
1957 3.650223
1958 3.127558
1959 3.339042
1960 3.637356
1961 2.694270
1962 3.846459
1963 3.872155
1964 3.495143
1965 3.563646
1966 3.484771
1967 3.535642
1968 3.149760
1969 3.324331
1970 3.910590
1971 4.181354
1972 3.643518
1973 3.221640
1974 3.429388
1975 3.262061
1976 3.252163
1977 3.211674
1978 3.143541
1979 3.811463
1980 3.391151
1981 3.484396
1982 3.300293
1983 3.280656
1984 3.926091
1985 3.341804
1986 3.646847
1987 3.023891
1988 3.368552
1989 3.656772
1990 3.124970
1991 3.395892
1992 3.601489
1993 3.184838
1994 3.336835
1995 3.620765
1996 3.591308
1997 3.665893
1998 2.970907
1999 3.527431
2000 3.178621
2001 3.175938
2002 3.217070
2003 3.307998
2004 3.782079
2005 3.386965
2006 3.443764
2007 3.897259
2008 4.018257
2009 3.171267
2010 3.339862
2011 3.726831
2012 3.531950
2013 4.461862
2014 3.532266
2015 3.271466
2016 3.184080
2017 3.437080
2018 2.970357
2019 3.253822
2020 3.654399
2021 3.375819
2022 3.482027 
ggplot(del_norte_standard_dev_all, aes(x = waterYear, y = sd_2))+#, color = waterYear)) +
  geom_line(size= 0.7) +
  #geom_line(aes) +
  theme_few() +
  geom_smooth(method = "lm", se=FALSE) +
  ylab('SD') + 
  xlab('Water year')

Full timeseries Mann-Kendall & Sen’s Slope

sd_mk_all <- mk.test(del_norte_standard_dev_all$sd_2)
print(sd_mk_all)
## 
##  Mann-Kendall trend test
## 
## data:  del_norte_standard_dev_all$sd_2
## z = 4.7174, n = 130, p-value = 2.389e-06
## alternative hypothesis: true S is not equal to 0
## sample estimates:
##            S         varS          tau 
## 2.345000e+03 2.468917e+05 2.796661e-01
sd_sens_all <- sens.slope(del_norte_standard_dev_all$sd_2)
print(sd_sens_all)
## 
##  Sen's slope
## 
## data:  del_norte_standard_dev_all$sd_2
## z = 4.7174, n = 130, p-value = 2.389e-06
## alternative hypothesis: true z is not equal to 0
## 95 percent confidence interval:
##  0.002192823 0.005139811
## sample estimates:
## Sen's slope 
##  0.00370004

Checking a random year from the timeseries.

del_norte_standard_dev_87 <- del_norte_standard_dev %>% 
  filter(waterYear == 1987) %>% 
  group_by(waterYear) %>% 
  mutate(resid_mean = mean(residual)) %>%
           mutate(sd_1 = residual-resid_mean)

del_norte_standard_dev_87 <- del_norte_standard_dev_87 %>%
  group_by(waterYear) %>%
  mutate(sd_2 = (((sum((sd_1)^2))/((sum(tabulate(del_norte_standard_dev_87$waterDay)))-1)))^(0.5))

head(del_norte_standard_dev_87$sd_2, 1)
## [1] 3.023891

Looks good.

Summer temperature standard deviation

del_norte_standard_dev_all_summer <- del_norte_standard_dev %>%
  filter(waterDay >= 244 & waterDay <= 335) %>% # this might be better off as daymonth rather than day of water year due to leap year
  group_by(waterYear) %>% 
  mutate(nmbr = n())

del_norte_standard_dev_all_summer <- del_norte_standard_dev_all_summer %>% 
  group_by(waterYear) %>% 
  mutate(resid_mean = mean(residual)) %>%
  mutate(sd_1 = residual-resid_mean) %>% 
  mutate(sd_2 = (((sum((sd_1)^2))/((nmbr-1))))^(0.5)) %>%
  distinct(sd_2, .keep_all = TRUE) %>% 
   select(waterYear, sd_2)

del_norte_standard_dev_all_summer %>% 
  kable(.,'html') %>%
  kable_styling() %>%
  scroll_box(width='250px',height='500px')
waterYear sd_2
1893 1.616702
1894 1.718305
1895 1.724864
1896 1.545493
1897 1.579815
1898 1.786310
1899 1.547891
1900 1.651407
1901 1.751337
1902 2.194414
1903 1.861475
1904 1.394076
1905 1.468259
1906 1.913958
1907 1.628179
1908 1.501797
1909 1.307557
1910 1.428952
1911 1.369587
1912 1.790366
1913 1.798927
1914 1.288604
1915 1.323734
1916 1.108404
1917 1.663379
1918 1.755367
1919 2.058061
1920 1.571214
1921 1.821923
1922 1.441676
1923 1.338064
1924 1.444628
1925 1.610093
1926 1.706964
1927 1.756002
1928 1.842276
1929 1.322699
1930 1.862230
1931 1.348963
1932 1.501092
1933 1.390344
1934 1.258517
1935 1.271380
1936 1.884833
1937 2.427488
1938 1.686255
1939 1.903632
1940 1.734099
1941 1.835093
1942 1.550206
1943 1.743478
1944 1.610272
1945 2.019325
1946 1.996323
1947 2.125673
1948 2.066701
1949 1.277370
1950 2.182231
1951 1.905283
1952 1.793921
1953 1.640122
1954 1.882011
1955 1.912099
1956 1.836740
1957 2.003684
1958 1.687804
1959 1.461210
1960 1.737903
1961 1.534011
1962 1.644398
1963 1.545998
1964 1.825389
1965 1.573278
1966 1.359337
1967 1.926409
1968 1.829141
1969 2.277293
1970 2.097444
1971 1.574299
1972 1.766186
1973 2.431246
1974 2.009254
1975 1.653547
1976 1.807136
1977 1.478783
1978 1.522808
1979 1.842975
1980 1.583710
1981 2.314608
1982 1.677585
1983 2.024668
1984 1.659605
1985 2.137010
1986 1.677283
1987 1.714587
1988 1.689609
1989 1.568497
1990 2.217696
1991 1.503807
1992 1.913593
1993 1.767483
1994 1.597623
1995 2.223239
1996 1.613703
1997 1.810866
1998 2.010522
1999 1.713251
2000 1.643233
2001 1.731585
2002 1.747722
2003 2.082463
2004 2.240902
2005 1.848249
2006 2.077114
2007 1.677673
2008 1.862177
2009 1.738860
2010 1.924970
2011 1.453084
2012 1.537989
2013 2.224635
2014 1.960436
2015 1.769018
2016 2.026886
2017 2.032470
2018 1.799639
2019 1.786262
2020 1.891842
2021 2.312794
2022 2.030081 
ggplot(del_norte_standard_dev_all_summer, aes(x = waterYear, y = sd_2))+#, color = waterYear)) +
  geom_line(size= 0.7) +
  #geom_line(aes) +
  theme_few() +
  geom_smooth(method = "lm", se=FALSE) +
  ylab('SD') + 
  xlab('Water year')

Jun-Aug standard deviation for water years 1900-2021

Mann-Kendall & Sen’s Slope

Summer standard deviations.

sd_mk_summer <- mk.test(del_norte_standard_dev_all_summer$sd_2)
print(sd_mk_summer)
## 
##  Mann-Kendall trend test
## 
## data:  del_norte_standard_dev_all_summer$sd_2
## z = 4.1499, n = 130, p-value = 3.327e-05
## alternative hypothesis: true S is not equal to 0
## sample estimates:
##            S         varS          tau 
## 2.063000e+03 2.468917e+05 2.460346e-01
sd_sens_summer <- sens.slope(del_norte_standard_dev_all_summer$sd_2)
print(sd_sens_summer)
## 
##  Sen's slope
## 
## data:  del_norte_standard_dev_all_summer$sd_2
## z = 4.1499, n = 130, p-value = 3.327e-05
## alternative hypothesis: true z is not equal to 0
## 95 percent confidence interval:
##  0.001370228 0.003733014
## sample estimates:
## Sen's slope 
## 0.002586023

Winter temperature standard deviation

del_norte_standard_dev_all_winter <- del_norte_standard_dev %>%
  filter(waterDay >= 32 & waterDay <= 182) %>% # this might be better off as daymonth rather than day of water year due to leap year
  group_by(waterYear) %>% 
  mutate(nmbr = n())

del_norte_standard_dev_all_winter <- del_norte_standard_dev_all_winter %>% 
  group_by(waterYear) %>% 
  mutate(resid_mean = mean(residual)) %>%
  mutate(sd_1 = residual-resid_mean) %>% 
  mutate(sd_2 = (((sum((sd_1)^2))/((nmbr-1))))^(0.5)) %>%
  distinct(sd_2, .keep_all = TRUE) %>% 
   select(waterYear, sd_2)

del_norte_standard_dev_all_winter %>% 
  kable(.,'html') %>%
  kable_styling() %>%
  scroll_box(width='250px',height='500px')
waterYear sd_2
1893 3.594948
1894 3.918018
1895 3.861162
1896 3.696494
1897 3.947023
1898 3.927829
1899 3.875649
1900 3.285578
1901 4.012814
1902 4.336653
1903 4.086321
1904 3.784659
1905 3.610353
1906 3.291979
1907 3.371937
1908 2.729394
1909 3.890538
1910 5.112050
1911 3.952553
1912 3.432775
1913 3.445698
1914 3.869463
1915 3.073329
1916 3.858286
1917 3.783862
1918 3.958981
1919 3.666477
1920 3.742279
1921 3.695989
1922 4.894954
1923 3.590838
1924 4.392766
1925 4.542657
1926 3.612205
1927 3.979158
1928 4.045089
1929 3.505447
1930 4.817291
1931 3.587287
1932 5.163963
1933 4.741613
1934 3.263905
1935 4.361538
1936 3.238379
1937 4.256182
1938 3.366052
1939 4.586916
1940 4.286941
1941 3.619731
1942 4.203736
1943 3.800891
1944 4.013060
1945 3.148753
1946 3.905203
1947 4.061710
1948 4.855627
1949 4.313098
1950 4.169221
1951 4.430057
1952 4.188744
1953 5.313154
1954 4.582321
1955 4.307055
1956 4.925458
1957 4.510967
1958 4.116403
1959 3.978053
1960 4.548065
1961 3.075771
1962 5.078738
1963 5.461302
1964 4.495260
1965 4.616120
1966 4.639273
1967 4.096578
1968 3.916700
1969 4.142396
1970 4.636417
1971 5.447052
1972 4.782661
1973 3.169461
1974 4.263881
1975 4.205329
1976 4.188874
1977 4.145999
1978 3.258374
1979 4.535424
1980 4.261592
1981 3.886644
1982 4.238960
1983 3.283920
1984 4.387854
1985 3.933623
1986 4.700489
1987 3.732172
1988 3.950937
1989 4.576652
1990 3.660709
1991 4.447670
1992 3.495688
1993 3.867566
1994 4.005916
1995 4.174198
1996 4.206461
1997 4.445972
1998 3.575336
1999 3.759359
2000 3.638909
2001 3.931654
2002 3.809549
2003 4.027876
2004 4.712365
2005 4.368520
2006 4.235769
2007 5.105147
2008 4.933385
2009 3.928699
2010 3.569531
2011 4.774349
2012 4.242636
2013 5.259241
2014 4.244943
2015 4.038394
2016 3.870281
2017 4.144996
2018 3.433675
2019 3.841134
2020 3.906713
2021 3.487736
2022 4.309173 
ggplot(del_norte_standard_dev_all_winter, aes(x = waterYear, y = sd_2))+#, color = waterYear)) +
  geom_line(size= 0.7) +
  #geom_line(aes) +
  theme_few() +
  geom_smooth(method = "lm", se=FALSE) +
  ylab('SD') + 
  xlab('Water year')

Nov-Mar standard deviation for water years 1900-2021

Mann-Kendall & Sen’s Slope

Winter standard deviations.

sd_mk_winter <- mk.test(del_norte_standard_dev_all_winter$sd_2)
print(sd_mk_winter)
## 
##  Mann-Kendall trend test
## 
## data:  del_norte_standard_dev_all_winter$sd_2
## z = 2.3305, n = 130, p-value = 0.01978
## alternative hypothesis: true S is not equal to 0
## sample estimates:
##            S         varS          tau 
## 1.159000e+03 2.468917e+05 1.382230e-01
sd_sens_winter <- sens.slope(del_norte_standard_dev_all_winter$sd_2)
print(sd_sens_winter)
## 
##  Sen's slope
## 
## data:  del_norte_standard_dev_all_winter$sd_2
## z = 2.3305, n = 130, p-value = 0.01978
## alternative hypothesis: true z is not equal to 0
## 95 percent confidence interval:
##  0.0005225851 0.0052967417
## sample estimates:
## Sen's slope 
## 0.002943756

Smaller time increments might be useful….

1987 - 2021

#average water year temperature

del_norte_yearly_wy_aver_87_21 <- del_norte_clean %>%
  filter(waterYear >= 1987 & waterYear <= 2021) %>% 
  group_by(waterYear) %>% 
  mutate(aver_ann_temp = mean(avg_T_c))

#Average temperature by day for all water years:

del_norte_daily_wy_aver_87_21 <- del_norte_yearly_wy_aver_87_21 %>% 
  group_by(daymonth) %>% 
  mutate(aver_day_temp = mean(avg_T_c))

#average mean temperature by day for the period of record:

del_norte_daily_wy_aver_87_21 <- del_norte_daily_wy_aver_87_21 %>% 
  group_by(daymonth) %>% 
  mutate(all_ave_temp = mean(del_norte_daily_wy_aver_87_21$aver_day_temp)) %>% 
  group_by(waterYear) %>% 
  mutate(nmbr = n())
  

del_norte_standard_dev_87_21 <- del_norte_daily_wy_aver_87_21 %>% 
  group_by(waterYear) %>% 
  mutate(residual = (all_ave_temp-aver_ann_temp)+avg_T_c-aver_day_temp) %>% 
  mutate(deviation = abs(residual-lag(residual)))

del_norte_standard_dev_87_21 <- del_norte_standard_dev_87_21 %>% 
  group_by(waterYear) %>% 
  mutate(resid_mean = mean(residual)) %>%
  mutate(sd_1 = residual-resid_mean) %>% 
  mutate(sd_2 = (((sum((sd_1)^2))/((nmbr-1))))^(0.5)) %>%
  distinct(sd_2, .keep_all = TRUE) %>% 
   select(waterYear, sd_2)


del_norte_standard_dev_87_21 %>% 
  kable(.,'html') %>%
  kable_styling() %>%
  scroll_box(width='250px',height='500px')
waterYear sd_2
1987 3.049446
1988 3.274043
1989 3.541523
1990 3.099801
1991 3.311386
1992 3.458707
1993 3.042800
1994 3.241010
1995 3.713253
1996 3.656668
1997 3.580422
1998 2.919295
1999 3.670739
2000 3.368763
2001 3.103862
2002 3.124365
2003 3.283259
2004 3.745359
2005 3.402699
2006 3.495618
2007 3.690314
2008 3.887032
2009 3.156753
2010 3.308721
2011 3.722723
2012 3.320809
2013 4.310065
2014 3.362813
2015 3.325711
2016 3.130585
2017 3.362125
2018 3.095342
2019 3.249459
2020 3.440448
2021 3.326127 
ggplot(del_norte_standard_dev_87_21, aes(x = waterYear, y = sd_2))+#, color = waterYear)) +
  geom_line(size= 0.7) +
  #geom_line(aes) +
  theme_few() +
  geom_smooth(method = "lm", se=FALSE) +
  ylab('SD') + 
  xlab('Water year')

SD for water years 1987 - 2021

1987 - 2021 Mann-Kendall & Sen’s Slope

sd_mk_87_21 <- mk.test(del_norte_standard_dev_87_21$sd_2)
print(sd_mk_87_21)
## 
##  Mann-Kendall trend test
## 
## data:  del_norte_standard_dev_87_21$sd_2
## z = 0.82368, n = 35, p-value = 0.4101
## alternative hypothesis: true S is not equal to 0
## sample estimates:
##            S         varS          tau 
## 5.900000e+01 4.958333e+03 9.915966e-02
sd_sens_87_21 <- sens.slope(del_norte_standard_dev_87_21$sd_2)
print(sd_sens_87_21)
## 
##  Sen's slope
## 
## data:  del_norte_standard_dev_87_21$sd_2
## z = 0.82368, n = 35, p-value = 0.4101
## alternative hypothesis: true z is not equal to 0
## 95 percent confidence interval:
##  -0.006279148  0.011354897
## sample estimates:
## Sen's slope 
## 0.003037291

Summer and Winter 87-21

Summer temperature standard deviation

# using the 1987- 2021 data frame

del_norte_standard_dev_all_87_21_summer <- del_norte_daily_wy_aver_87_21 %>%
  filter(waterDay >= 244 & waterDay <= 335) %>% # this might be better off as daymonth rather than day of water year due to leap year
  group_by(waterYear) %>% 
  mutate(residual = (all_ave_temp-aver_ann_temp)+avg_T_c-aver_day_temp) %>% 
  mutate(deviation = abs(residual-lag(residual))) %>% 
  group_by(waterYear) %>% 
  mutate(nmbr = n())

del_norte_standard_dev_all_87_21_summer <- del_norte_standard_dev_all_87_21_summer %>% 
  group_by(waterYear) %>% 
  mutate(resid_mean = mean(residual)) %>%
  mutate(sd_1 = residual-resid_mean) %>% 
  mutate(sd_2 = (((sum((sd_1)^2))/((nmbr-1))))^(0.5)) %>%
  distinct(sd_2, .keep_all = TRUE) %>% 
   select(waterYear, sd_2)

del_norte_standard_dev_all_87_21_summer %>% 
  kable(.,'html') %>%
  kable_styling() %>%
  scroll_box(width='250px',height='500px')
waterYear sd_2
1987 1.666876
1988 1.638926
1989 1.565457
1990 2.124272
1991 1.542094
1992 1.881595
1993 1.830733
1994 1.560719
1995 2.318927
1996 1.542480
1997 1.722056
1998 1.983380
1999 1.738525
2000 1.628533
2001 1.634000
2002 1.675954
2003 2.113700
2004 2.212596
2005 1.832730
2006 1.958719
2007 1.710145
2008 1.857423
2009 1.784096
2010 1.871881
2011 1.483893
2012 1.475483
2013 2.210304
2014 1.848449
2015 1.745719
2016 1.954333
2017 1.915606
2018 1.770882
2019 1.863582
2020 1.875031
2021 2.348222 
ggplot(del_norte_standard_dev_all_87_21_summer, aes(x = waterYear, y = sd_2))+#, color = waterYear)) +
  geom_line(size= 0.7) +
  #geom_line(aes) +
  theme_few() +
  geom_smooth(method = "lm", se=FALSE) +
  ylab('SD') + 
  xlab('Water year')

Summer season SD for water years 1987 - 2021

Mann-Kendall & Sen’s Slope

Summer 87-21 standard deviations.

sd_mk_summer <- mk.test(del_norte_standard_dev_all_87_21_summer$sd_2)
print(sd_mk_summer)
## 
##  Mann-Kendall trend test
## 
## data:  del_norte_standard_dev_all_87_21_summer$sd_2
## z = 1.7326, n = 35, p-value = 0.08317
## alternative hypothesis: true S is not equal to 0
## sample estimates:
##            S         varS          tau 
##  123.0000000 4958.3333333    0.2067227
sd_sens_summer <- sens.slope(del_norte_standard_dev_all_summer$sd_2)
print(sd_sens_summer)
## 
##  Sen's slope
## 
## data:  del_norte_standard_dev_all_summer$sd_2
## z = 4.1499, n = 130, p-value = 3.327e-05
## alternative hypothesis: true z is not equal to 0
## 95 percent confidence interval:
##  0.001370228 0.003733014
## sample estimates:
## Sen's slope 
## 0.002586023

Winter temperature standard deviation

# using the 1987- 2021 data frame

del_norte_standard_dev_all_87_21_winter <- del_norte_daily_wy_aver_87_21 %>%
  filter(waterDay >= 32 & waterDay <= 182) %>% # this might be better off as daymonth rather than day of water year due to leap year
  group_by(waterYear) %>% 
  mutate(residual = (all_ave_temp-aver_ann_temp)+avg_T_c-aver_day_temp) %>% 
  mutate(deviation = abs(residual-lag(residual))) %>% 
  group_by(waterYear) %>% 
  mutate(nmbr = n())

del_norte_standard_dev_all_87_21_winter <- del_norte_standard_dev_all_87_21_winter %>% 
  group_by(waterYear) %>% 
  mutate(resid_mean = mean(residual)) %>%
  mutate(sd_1 = residual-resid_mean) %>% 
  mutate(sd_2 = (((sum((sd_1)^2))/((nmbr-1))))^(0.5)) %>%
  distinct(sd_2, .keep_all = TRUE) %>% 
   select(waterYear, sd_2)

del_norte_standard_dev_all_87_21_winter %>% 
  kable(.,'html') %>%
  kable_styling() %>%
  scroll_box(width='250px',height='500px')
waterYear sd_2
1987 3.813461
1988 3.892916
1989 4.404551
1990 3.573400
1991 4.316134
1992 3.366099
1993 3.640254
1994 3.883465
1995 4.179941
1996 4.277382
1997 4.256606
1998 3.481668
1999 3.837016
2000 3.877807
2001 3.904396
2002 3.811252
2003 3.948204
2004 4.612768
2005 4.335625
2006 4.355642
2007 4.797987
2008 4.829273
2009 3.841686
2010 3.665058
2011 4.831621
2012 4.009317
2013 5.176774
2014 4.002662
2015 4.041898
2016 3.764093
2017 3.923724
2018 3.631138
2019 3.914099
2020 3.678198
2021 3.514880 
ggplot(del_norte_standard_dev_all_87_21_winter, aes(x = waterYear, y = sd_2))+#, color = waterYear)) +
  geom_line(size= 0.7) +
  #geom_line(aes) +
  theme_few() +
  geom_smooth(method = "lm", se=FALSE) +
  ylab('SD') + 
  xlab('Water year')

Winter season SD for water years 1987 - 2021

Mann-Kendall & Sen’s Slope

1987 - 2021 winter Mann-Kendall & Sen’s Slope

sd_mk_87_21_winter <- mk.test(del_norte_standard_dev_all_87_21_winter$sd_2)
print(sd_mk_87_21_winter)
## 
##  Mann-Kendall trend test
## 
## data:  del_norte_standard_dev_all_87_21_winter$sd_2
## z = 0.31243, n = 35, p-value = 0.7547
## alternative hypothesis: true S is not equal to 0
## sample estimates:
##            S         varS          tau 
## 2.300000e+01 4.958333e+03 3.865546e-02
sd_sens_87_21_winter <- sens.slope(del_norte_standard_dev_all_87_21_winter$sd_2)
print(sd_sens_87_21_winter)
## 
##  Sen's slope
## 
## data:  del_norte_standard_dev_all_87_21_winter$sd_2
## z = 0.31243, n = 35, p-value = 0.7547
## alternative hypothesis: true z is not equal to 0
## 95 percent confidence interval:
##  -0.01360253  0.01569025
## sample estimates:
## Sen's slope 
## 0.001405344

1900-1915 minimum and maximum temperatures

early_del_norte <- del_norte_clean %>% 
  filter(waterYear >= 1893 & waterYear <= 1915)

del_norte_temp_min_xts <- xts(early_del_norte$lowc, order.by = early_del_norte$Date)

dygraph(del_norte_temp_min_xts) %>%
  dyAxis("y", label = "Daily min temperature (°C)")

Daily minimum temperatures for water years 1900-1915

del_norte_temp_max_xts <- xts(early_del_norte$highc, order.by = early_del_norte$Date)

dygraph(del_norte_temp_max_xts) %>%
  dyAxis("y", label = "Daily max temperature (°C)")

Daily minimum temperatures for water years 1900-1915

By groups of decades

1900-1930

#average water year temperature

del_norte_yearly_wy_aver_00_30 <- del_norte_clean %>%
  filter(waterYear >= 1900 & waterYear <= 1930) %>% 
  group_by(waterYear) %>% 
  mutate(aver_ann_temp = mean(avg_T_c))

#Average temperature by day for all water years:

del_norte_daily_wy_aver_00_30 <- del_norte_yearly_wy_aver_00_30 %>% 
  group_by(daymonth) %>% 
  mutate(aver_day_temp = mean(avg_T_c))

#average mean temperature by day for the period of record:

del_norte_daily_wy_aver_00_30 <- del_norte_daily_wy_aver_00_30 %>% 
  group_by(daymonth) %>% 
  mutate(all_ave_temp = mean(del_norte_daily_wy_aver_00_30$aver_day_temp)) %>% 
  group_by(waterYear) %>% 
  mutate(nmbr = n())
  

del_norte_standard_dev_00_30 <- del_norte_daily_wy_aver_00_30 %>% 
  group_by(waterYear) %>% 
  mutate(residual = (all_ave_temp-aver_ann_temp)+avg_T_c-aver_day_temp) %>% 
  mutate(deviation = abs(residual-lag(residual)))

del_norte_standard_dev_00_30 <- del_norte_standard_dev_00_30 %>% 
  group_by(waterYear) %>% 
  mutate(resid_mean = mean(residual)) %>%
  mutate(sd_1 = residual-resid_mean) %>% 
  mutate(sd_2 = (((sum((sd_1)^2))/((nmbr-1))))^(0.5)) %>%
  distinct(sd_2, .keep_all = TRUE) %>% 
   select(waterYear, sd_2)


del_norte_standard_dev_00_30 %>% 
  kable(.,'html') %>%
  kable_styling() %>%
  scroll_box(width='250px',height='500px')
waterYear sd_2
1900 2.670640
1901 3.048471
1902 3.170437
1903 3.265402
1904 2.847399
1905 2.714680
1906 2.652614
1907 3.505271
1908 2.402912
1909 3.184523
1910 3.698099
1911 3.166985
1912 2.802507
1913 2.970735
1914 2.921311
1915 2.583496
1916 2.912074
1917 3.079872
1918 3.150513
1919 3.739201
1920 3.145837
1921 3.072360
1922 3.580464
1923 2.578237
1924 3.233255
1925 3.417131
1926 2.894750
1927 3.398953
1928 3.422627
1929 2.791312
1930 3.644880 
ggplot(del_norte_standard_dev_00_30, aes(x = waterYear, y = sd_2))+#, color = waterYear)) +
  geom_line(size= 0.7) +
  #geom_line(aes) +
  theme_few() +
  geom_smooth(method = "lm", se=FALSE) +
  ylab('SD') + 
  xlab('Water year')

SD for water years 1900 - 1930

1900 - 1930 Mann-Kendall & Sen’s Slope

sd_mk_00_30 <- mk.test(del_norte_standard_dev_00_30$sd_2)
print(sd_mk_00_30)
## 
##  Mann-Kendall trend test
## 
## data:  del_norte_standard_dev_00_30$sd_2
## z = 1.4617, n = 31, p-value = 0.1438
## alternative hypothesis: true S is not equal to 0
## sample estimates:
##            S         varS          tau 
##   87.0000000 3461.6666667    0.1870968
sd_sens_00_30 <- sens.slope(del_norte_standard_dev_00_30$sd_2)
print(sd_sens_00_30)
## 
##  Sen's slope
## 
## data:  del_norte_standard_dev_00_30$sd_2
## z = 1.4617, n = 31, p-value = 0.1438
## alternative hypothesis: true z is not equal to 0
## 95 percent confidence interval:
##  -0.003004377  0.029983672
## sample estimates:
## Sen's slope 
##  0.01405475

Summer and Winter 00-30

Summer temperature standard deviation

# using the 1900- 1930 data frame

del_norte_standard_dev_all_00_30_summer <- del_norte_daily_wy_aver_00_30 %>%
  filter(waterDay >= 244 & waterDay <= 335) %>% # this might be better off as daymonth rather than day of water year due to leap year
  group_by(waterYear) %>% 
  mutate(residual = (all_ave_temp-aver_ann_temp)+avg_T_c-aver_day_temp) %>% 
  mutate(deviation = abs(residual-lag(residual))) %>% 
  group_by(waterYear) %>% 
  mutate(nmbr = n())

del_norte_standard_dev_all_00_30_summer <- del_norte_standard_dev_all_00_30_summer %>% 
  group_by(waterYear) %>% 
  mutate(resid_mean = mean(residual)) %>%
  mutate(sd_1 = residual-resid_mean) %>% 
  mutate(sd_2 = (((sum((sd_1)^2))/((nmbr-1))))^(0.5)) %>%
  distinct(sd_2, .keep_all = TRUE) %>% 
   select(waterYear, sd_2)

del_norte_standard_dev_all_00_30_summer %>% 
  kable(.,'html') %>%
  kable_styling() %>%
  scroll_box(width='250px',height='500px')
waterYear sd_2
1900 1.658110
1901 1.649685
1902 2.128479
1903 1.788408
1904 1.387950
1905 1.465325
1906 1.864438
1907 1.646117
1908 1.458500
1909 1.286226
1910 1.440182
1911 1.461181
1912 1.725176
1913 1.857125
1914 1.280254
1915 1.253525
1916 1.104074
1917 1.576149
1918 1.792958
1919 2.011147
1920 1.567232
1921 1.714232
1922 1.435451
1923 1.282250
1924 1.437670
1925 1.558668
1926 1.687437
1927 1.768045
1928 1.856520
1929 1.251370
1930 1.919901 
ggplot(del_norte_standard_dev_all_00_30_summer, aes(x = waterYear, y = sd_2))+#, color = waterYear)) +
  geom_line(size= 0.7) +
  #geom_line(aes) +
  theme_few() +
  geom_smooth(method = "lm", se=FALSE) +
  ylab('SD') + 
  xlab('Water year')

Summer season SD for water years 1900 - 1930

Mann-Kendall & Sen’s Slope

Summer 00-30 standard deviations.

sd_mk_00_30_summer <- mk.test(del_norte_standard_dev_all_00_30_summer$sd_2)
print(sd_mk_00_30_summer)
## 
##  Mann-Kendall trend test
## 
## data:  del_norte_standard_dev_all_00_30_summer$sd_2
## z = -0.30594, n = 31, p-value = 0.7597
## alternative hypothesis: true S is not equal to 0
## sample estimates:
##             S          varS           tau 
##  -19.00000000 3461.66666667   -0.04086022
sd_sens_00_30_summer <- sens.slope(del_norte_standard_dev_all_00_30_summer$sd_2)
print(sd_sens_00_30_summer)
## 
##  Sen's slope
## 
## data:  del_norte_standard_dev_all_00_30_summer$sd_2
## z = -0.30594, n = 31, p-value = 0.7597
## alternative hypothesis: true z is not equal to 0
## 95 percent confidence interval:
##  -0.01175000  0.01120514
## sample estimates:
##  Sen's slope 
## -0.001044823

Winter temperature standard deviation

# using the 00-30 data frame

del_norte_standard_dev_all_00_30_winter <- del_norte_daily_wy_aver_00_30 %>%
  filter(waterDay >= 32 & waterDay <= 182) %>% # this might be better off as daymonth rather than day of water year due to leap year
  group_by(waterYear) %>% 
  mutate(residual = (all_ave_temp-aver_ann_temp)+avg_T_c-aver_day_temp) %>% 
  mutate(deviation = abs(residual-lag(residual))) %>% 
  group_by(waterYear) %>% 
  mutate(nmbr = n())

del_norte_standard_dev_all_00_30_winter <- del_norte_standard_dev_all_00_30_winter %>% 
  group_by(waterYear) %>% 
  mutate(resid_mean = mean(residual)) %>%
  mutate(sd_1 = residual-resid_mean) %>% 
  mutate(sd_2 = (((sum((sd_1)^2))/((nmbr-1))))^(0.5)) %>%
  distinct(sd_2, .keep_all = TRUE) %>% 
   select(waterYear, sd_2)

del_norte_standard_dev_all_00_30_winter %>% 
  kable(.,'html') %>%
  kable_styling() %>%
  scroll_box(width='250px',height='500px')
waterYear sd_2
1900 3.179435
1901 3.915741
1902 4.137566
1903 4.119047
1904 3.722432
1905 3.586544
1906 3.310207
1907 3.461982
1908 2.748450
1909 3.731511
1910 5.013542
1911 4.025188
1912 3.433575
1913 3.314971
1914 3.782262
1915 3.148231
1916 3.759962
1917 3.659457
1918 3.953960
1919 3.891303
1920 3.634187
1921 3.610989
1922 4.945866
1923 3.512318
1924 4.340311
1925 4.485393
1926 3.585932
1927 3.728669
1928 3.834410
1929 3.473418
1930 4.760612 
ggplot(del_norte_standard_dev_all_00_30_winter, aes(x = waterYear, y = sd_2))+#, color = waterYear)) +
  geom_line(size= 0.7) +
  #geom_line(aes) +
  theme_few() +
  geom_smooth(method = "lm", se=FALSE) +
  ylab('SD') + 
  xlab('Water year')

Winter season SD for water years 1900 - 1930

Mann-Kendall & Sen’s Slope

1900 - 1930 winter Mann-Kendall & Sen’s Slope

sd_mk_00_30_winter <- mk.test(del_norte_standard_dev_all_00_30_winter$sd_2)
print(sd_mk_00_30_winter)
## 
##  Mann-Kendall trend test
## 
## data:  del_norte_standard_dev_all_00_30_winter$sd_2
## z = 0.88381, n = 31, p-value = 0.3768
## alternative hypothesis: true S is not equal to 0
## sample estimates:
##            S         varS          tau 
##   53.0000000 3461.6666667    0.1139785
sd_sens_00_30_winter <- sens.slope(del_norte_standard_dev_all_00_30_winter$sd_2)
print(sd_sens_00_30_winter)
## 
##  Sen's slope
## 
## data:  del_norte_standard_dev_all_00_30_winter$sd_2
## z = 0.88381, n = 31, p-value = 0.3768
## alternative hypothesis: true z is not equal to 0
## 95 percent confidence interval:
##  -0.01105860  0.03191729
## sample estimates:
## Sen's slope 
##  0.01042583

1931-1960

#average water year temperature

del_norte_yearly_wy_aver_31_60 <- del_norte_clean %>%
  filter(waterYear >= 1931 & waterYear <= 1960) %>% 
  group_by(waterYear) %>% 
  mutate(aver_ann_temp = mean(avg_T_c))

#Average temperature by day for all water years:

del_norte_daily_wy_aver_31_60 <- del_norte_yearly_wy_aver_31_60 %>% 
  group_by(daymonth) %>% 
  mutate(aver_day_temp = mean(avg_T_c))

#average mean temperature by day for the period of record:

del_norte_daily_wy_aver_31_60 <- del_norte_daily_wy_aver_31_60 %>% 
  group_by(daymonth) %>% 
  mutate(all_ave_temp = mean(del_norte_daily_wy_aver_31_60$aver_day_temp)) %>% 
  group_by(waterYear) %>% 
  mutate(nmbr = n())
  

del_norte_standard_dev_31_60 <- del_norte_daily_wy_aver_31_60 %>% 
  group_by(waterYear) %>% 
  mutate(residual = (all_ave_temp-aver_ann_temp)+avg_T_c-aver_day_temp) %>% 
  mutate(deviation = abs(residual-lag(residual)))

del_norte_standard_dev_31_60 <- del_norte_standard_dev_31_60 %>% 
  group_by(waterYear) %>% 
  mutate(resid_mean = mean(residual)) %>%
  mutate(sd_1 = residual-resid_mean) %>% 
  mutate(sd_2 = (((sum((sd_1)^2))/((nmbr-1))))^(0.5)) %>%
  distinct(sd_2, .keep_all = TRUE) %>% 
   select(waterYear, sd_2)


del_norte_standard_dev_31_60 %>% 
  kable(.,'html') %>%
  kable_styling() %>%
  scroll_box(width='250px',height='500px')
waterYear sd_2
1931 2.885049
1932 3.770130
1933 3.852008
1934 2.838981
1935 3.380770
1936 2.872823
1937 3.309781
1938 3.071443
1939 3.587949
1940 3.077642
1941 3.191093
1942 3.047844
1943 3.060907
1944 3.193958
1945 2.851190
1946 3.120056
1947 3.301030
1948 3.802148
1949 3.270271
1950 3.504621
1951 3.396156
1952 3.288819
1953 4.020477
1954 3.389217
1955 3.279470
1956 3.514359
1957 3.582358
1958 3.051335
1959 3.255626
1960 3.541445 
ggplot(del_norte_standard_dev_31_60, aes(x = waterYear, y = sd_2))+#, color = waterYear)) +
  geom_line(size= 0.7) +
  #geom_line(aes) +
  theme_few() +
  geom_smooth(method = "lm", se=FALSE) +
  ylab('SD') + 
  xlab('Water year')

SD for water years 1931 - 1960

1931 - 1960 Mann-Kendall & Sen’s Slope

sd_mk_31_60 <- mk.test(del_norte_standard_dev_31_60$sd_2)
print(sd_mk_31_60)
## 
##  Mann-Kendall trend test
## 
## data:  del_norte_standard_dev_31_60$sd_2
## z = 1.57, n = 30, p-value = 0.1164
## alternative hypothesis: true S is not equal to 0
## sample estimates:
##            S         varS          tau 
##   89.0000000 3141.6666667    0.2045977
sd_sens_31_60 <- sens.slope(del_norte_standard_dev_31_60$sd_2)
print(sd_sens_31_60)
## 
##  Sen's slope
## 
## data:  del_norte_standard_dev_31_60$sd_2
## z = 1.57, n = 30, p-value = 0.1164
## alternative hypothesis: true z is not equal to 0
## 95 percent confidence interval:
##  -0.002312981  0.024977931
## sample estimates:
## Sen's slope 
##  0.01345518

Summer and Winter 31-60

Summer temperature standard deviation

# using the 1931- 1960 data frame

del_norte_standard_dev_all_31_60_summer <- del_norte_daily_wy_aver_31_60 %>%
  filter(waterDay >= 244 & waterDay <= 335) %>% # this might be better off as daymonth rather than day of water year due to leap year
  group_by(waterYear) %>% 
  mutate(residual = (all_ave_temp-aver_ann_temp)+avg_T_c-aver_day_temp) %>% 
  mutate(deviation = abs(residual-lag(residual))) %>% 
  group_by(waterYear) %>% 
  mutate(nmbr = n())

del_norte_standard_dev_all_31_60_summer <- del_norte_standard_dev_all_31_60_summer %>% 
  group_by(waterYear) %>% 
  mutate(resid_mean = mean(residual)) %>%
  mutate(sd_1 = residual-resid_mean) %>% 
  mutate(sd_2 = (((sum((sd_1)^2))/((nmbr-1))))^(0.5)) %>%
  distinct(sd_2, .keep_all = TRUE) %>% 
   select(waterYear, sd_2)

del_norte_standard_dev_all_31_60_summer %>% 
  kable(.,'html') %>%
  kable_styling() %>%
  scroll_box(width='250px',height='500px')
waterYear sd_2
1931 1.357917
1932 1.444201
1933 1.441354
1934 1.272179
1935 1.209747
1936 1.837705
1937 2.354713
1938 1.630502
1939 1.946711
1940 1.652860
1941 1.743069
1942 1.523017
1943 1.703816
1944 1.519863
1945 1.988650
1946 2.014436
1947 2.096767
1948 2.069391
1949 1.283852
1950 2.173349
1951 1.866589
1952 1.773101
1953 1.646071
1954 1.856894
1955 1.877769
1956 1.840758
1957 2.065215
1958 1.685443
1959 1.456781
1960 1.693242 
ggplot(del_norte_standard_dev_all_31_60_summer, aes(x = waterYear, y = sd_2))+#, color = waterYear)) +
  geom_line(size= 0.7) +
  #geom_line(aes) +
  theme_few() +
  geom_smooth(method = "lm", se=FALSE) +
  ylab('SD') + 
  xlab('Water year')

Summer season SD for water years 1931 - 1960

Mann-Kendall & Sen’s Slope

Summer 31-60 standard deviations.

sd_mk_31_60_summer <- mk.test(del_norte_standard_dev_all_31_60_summer$sd_2)
print(sd_mk_31_60_summer)
## 
##  Mann-Kendall trend test
## 
## data:  del_norte_standard_dev_all_31_60_summer$sd_2
## z = 1.7127, n = 30, p-value = 0.08676
## alternative hypothesis: true S is not equal to 0
## sample estimates:
##            S         varS          tau 
##   97.0000000 3141.6666667    0.2229885
sd_sens_31_60_summer <- sens.slope(del_norte_standard_dev_all_31_60_summer$sd_2)
print(sd_sens_31_60_summer)
## 
##  Sen's slope
## 
## data:  del_norte_standard_dev_all_31_60_summer$sd_2
## z = 1.7127, n = 30, p-value = 0.08676
## alternative hypothesis: true z is not equal to 0
## 95 percent confidence interval:
##  -0.002846981  0.023974556
## sample estimates:
## Sen's slope 
##  0.01118667

Winter temperature standard deviation

# using the 31-60 data frame

del_norte_standard_dev_all_31_60_winter <- del_norte_daily_wy_aver_31_60 %>%
  filter(waterDay >= 32 & waterDay <= 182) %>% # this might be better off as daymonth rather than day of water year due to leap year
  group_by(waterYear) %>% 
  mutate(residual = (all_ave_temp-aver_ann_temp)+avg_T_c-aver_day_temp) %>% 
  mutate(deviation = abs(residual-lag(residual))) %>% 
  group_by(waterYear) %>% 
  mutate(nmbr = n())

del_norte_standard_dev_all_31_60_winter <- del_norte_standard_dev_all_31_60_winter %>% 
  group_by(waterYear) %>% 
  mutate(resid_mean = mean(residual)) %>%
  mutate(sd_1 = residual-resid_mean) %>% 
  mutate(sd_2 = (((sum((sd_1)^2))/((nmbr-1))))^(0.5)) %>%
  distinct(sd_2, .keep_all = TRUE) %>% 
   select(waterYear, sd_2)

del_norte_standard_dev_all_31_60_winter %>% 
  kable(.,'html') %>%
  kable_styling() %>%
  scroll_box(width='250px',height='500px')
waterYear sd_2
1931 3.644805
1932 5.170104
1933 4.947210
1934 3.322783
1935 4.275249
1936 3.149329
1937 4.163917
1938 3.449850
1939 4.597173
1940 4.270947
1941 3.683929
1942 4.016393
1943 3.641176
1944 4.104933
1945 3.017130
1946 3.902809
1947 3.883944
1948 4.671679
1949 4.237985
1950 4.277180
1951 4.347027
1952 4.047289
1953 5.283903
1954 4.563898
1955 4.265180
1956 4.730117
1957 4.391196
1958 3.990812
1959 3.813852
1960 4.369029 
ggplot(del_norte_standard_dev_all_31_60_winter, aes(x = waterYear, y = sd_2))+#, color = waterYear)) +
  geom_line(size= 0.7) +
  #geom_line(aes) +
  theme_few() +
  geom_smooth(method = "lm", se=FALSE) +
  ylab('SD') + 
  xlab('Water year')

Winter season SD for water years 1931 - 1960

Mann-Kendall & Sen’s Slope

1931 - 1960 winter Mann-Kendall & Sen’s Slope

sd_mk_31_60_winter <- mk.test(del_norte_standard_dev_all_31_60_winter$sd_2)
print(sd_mk_31_60_winter)
## 
##  Mann-Kendall trend test
## 
## data:  del_norte_standard_dev_all_31_60_winter$sd_2
## z = 1.2489, n = 30, p-value = 0.2117
## alternative hypothesis: true S is not equal to 0
## sample estimates:
##            S         varS          tau 
##   71.0000000 3141.6666667    0.1632184
sd_sens_31_60_winter <- sens.slope(del_norte_standard_dev_all_31_60_winter$sd_2)
print(sd_sens_31_60_winter)
## 
##  Sen's slope
## 
## data:  del_norte_standard_dev_all_31_60_winter$sd_2
## z = 1.2489, n = 30, p-value = 0.2117
## alternative hypothesis: true z is not equal to 0
## 95 percent confidence interval:
##  -0.01086401  0.04151798
## sample estimates:
## Sen's slope 
##  0.01650605

1961-1990

#average water year temperature

del_norte_yearly_wy_aver_61_90 <- del_norte_clean %>%
  filter(waterYear >= 1961 & waterYear <= 1990) %>% 
  group_by(waterYear) %>% 
  mutate(aver_ann_temp = mean(avg_T_c))

#Average temperature by day for all water years:

del_norte_daily_wy_aver_61_90 <- del_norte_yearly_wy_aver_61_90 %>% 
  group_by(daymonth) %>% 
  mutate(aver_day_temp = mean(avg_T_c))

#average mean temperature by day for the period of record:

del_norte_daily_wy_aver_61_90 <- del_norte_daily_wy_aver_61_90 %>% 
  group_by(daymonth) %>% 
  mutate(all_ave_temp = mean(del_norte_daily_wy_aver_61_90$aver_day_temp)) %>% 
  group_by(waterYear) %>% 
  mutate(nmbr = n())
  

del_norte_standard_dev_61_90 <- del_norte_daily_wy_aver_61_90 %>% 
  group_by(waterYear) %>% 
  mutate(residual = (all_ave_temp-aver_ann_temp)+avg_T_c-aver_day_temp) %>% 
  mutate(deviation = abs(residual-lag(residual)))

del_norte_standard_dev_61_90 <- del_norte_standard_dev_61_90 %>% 
  group_by(waterYear) %>% 
  mutate(resid_mean = mean(residual)) %>%
  mutate(sd_1 = residual-resid_mean) %>% 
  mutate(sd_2 = (((sum((sd_1)^2))/((nmbr-1))))^(0.5)) %>%
  distinct(sd_2, .keep_all = TRUE) %>% 
   select(waterYear, sd_2)


del_norte_standard_dev_61_90 %>% 
  kable(.,'html') %>%
  kable_styling() %>%
  scroll_box(width='250px',height='500px')
waterYear sd_2
1961 2.652011
1962 3.838105
1963 3.849170
1964 3.400975
1965 3.510945
1966 3.457497
1967 3.418993
1968 3.074156
1969 3.359801
1970 3.825266
1971 4.060860
1972 3.569781
1973 3.170185
1974 3.312662
1975 3.157682
1976 3.102911
1977 3.190448
1978 3.131516
1979 3.771941
1980 3.361648
1981 3.511685
1982 3.252683
1983 3.364953
1984 3.869026
1985 3.356020
1986 3.555517
1987 3.118324
1988 3.251738
1989 3.520066
1990 3.140466 
ggplot(del_norte_standard_dev_61_90, aes(x = waterYear, y = sd_2))+#, color = waterYear)) +
  geom_line(size= 0.7) +
  #geom_line(aes) +
  theme_few() +
  geom_smooth(method = "lm", se=FALSE) +
  ylab('SD') + 
  xlab('Water year')

SD for water years 1961 - 1990

1961 - 1990 Mann-Kendall & Sen’s Slope

sd_mk_61_90 <- mk.test(del_norte_standard_dev_61_90$sd_2)
print(sd_mk_61_90)
## 
##  Mann-Kendall trend test
## 
## data:  del_norte_standard_dev_61_90$sd_2
## z = -0.82069, n = 30, p-value = 0.4118
## alternative hypothesis: true S is not equal to 0
## sample estimates:
##           S        varS         tau 
##  -47.000000 3141.666667   -0.108046
sd_sens_61_90 <- sens.slope(del_norte_standard_dev_61_90$sd_2)
print(sd_sens_61_90)
## 
##  Sen's slope
## 
## data:  del_norte_standard_dev_61_90$sd_2
## z = -0.82069, n = 30, p-value = 0.4118
## alternative hypothesis: true z is not equal to 0
## 95 percent confidence interval:
##  -0.019247118  0.007185514
## sample estimates:
##  Sen's slope 
## -0.005340888

Summer and Winter 61-90

Summer temperature standard deviation

# using the 1961- 1990 data frame

del_norte_standard_dev_all_61_90_summer <- del_norte_daily_wy_aver_61_90 %>%
  filter(waterDay >= 244 & waterDay <= 335) %>% # this might be better off as daymonth rather than day of water year due to leap year
  group_by(waterYear) %>% 
  mutate(residual = (all_ave_temp-aver_ann_temp)+avg_T_c-aver_day_temp) %>% 
  mutate(deviation = abs(residual-lag(residual))) %>% 
  group_by(waterYear) %>% 
  mutate(nmbr = n())

del_norte_standard_dev_all_61_90_summer <- del_norte_standard_dev_all_61_90_summer %>% 
  group_by(waterYear) %>% 
  mutate(resid_mean = mean(residual)) %>%
  mutate(sd_1 = residual-resid_mean) %>% 
  mutate(sd_2 = (((sum((sd_1)^2))/((nmbr-1))))^(0.5)) %>%
  distinct(sd_2, .keep_all = TRUE) %>% 
   select(waterYear, sd_2)

del_norte_standard_dev_all_61_90_summer %>% 
  kable(.,'html') %>%
  kable_styling() %>%
  scroll_box(width='250px',height='500px')
waterYear sd_2
1961 1.552531
1962 1.563890
1963 1.471129
1964 1.781435
1965 1.562627
1966 1.327834
1967 1.869669
1968 1.829171
1969 2.254014
1970 2.047500
1971 1.573143
1972 1.842777
1973 2.353221
1974 1.969620
1975 1.575966
1976 1.778709
1977 1.501636
1978 1.556862
1979 1.782639
1980 1.538799
1981 2.258635
1982 1.740249
1983 2.010755
1984 1.651535
1985 2.118872
1986 1.652165
1987 1.792312
1988 1.704304
1989 1.555450
1990 2.128402 
ggplot(del_norte_standard_dev_all_61_90_summer, aes(x = waterYear, y = sd_2))+#, color = waterYear)) +
  geom_line(size= 0.7) +
  #geom_line(aes) +
  theme_few() +
  geom_smooth(method = "lm", se=FALSE) +
  ylab('SD') + 
  xlab('Water year')

Summer season SD for water years 1961 - 1990

Mann-Kendall & Sen’s Slope

Summer 61-90 standard deviations.

sd_mk_61_90_summer <- mk.test(del_norte_standard_dev_all_61_90_summer$sd_2)
print(sd_mk_61_90_summer)
## 
##  Mann-Kendall trend test
## 
## data:  del_norte_standard_dev_all_61_90_summer$sd_2
## z = 1.1418, n = 30, p-value = 0.2535
## alternative hypothesis: true S is not equal to 0
## sample estimates:
##            S         varS          tau 
##   65.0000000 3141.6666667    0.1494253
sd_sens_61_90_summer <- sens.slope(del_norte_standard_dev_all_61_90_summer$sd_2)
print(sd_sens_61_90_summer)
## 
##  Sen's slope
## 
## data:  del_norte_standard_dev_all_61_90_summer$sd_2
## z = 1.1418, n = 30, p-value = 0.2535
## alternative hypothesis: true z is not equal to 0
## 95 percent confidence interval:
##  -0.004230137  0.016068404
## sample estimates:
## Sen's slope 
## 0.004758143

Winter temperature standard deviation

# using the 61-90 data frame

del_norte_standard_dev_all_61_90_winter <- del_norte_daily_wy_aver_61_90 %>%
  filter(waterDay >= 32 & waterDay <= 182) %>% # this might be better off as daymonth rather than day of water year due to leap year
  group_by(waterYear) %>% 
  mutate(residual = (all_ave_temp-aver_ann_temp)+avg_T_c-aver_day_temp) %>% 
  mutate(deviation = abs(residual-lag(residual))) %>% 
  group_by(waterYear) %>% 
  mutate(nmbr = n())

del_norte_standard_dev_all_61_90_winter <- del_norte_standard_dev_all_61_90_winter %>% 
  group_by(waterYear) %>% 
  mutate(resid_mean = mean(residual)) %>%
  mutate(sd_1 = residual-resid_mean) %>% 
  mutate(sd_2 = (((sum((sd_1)^2))/((nmbr-1))))^(0.5)) %>%
  distinct(sd_2, .keep_all = TRUE) %>% 
   select(waterYear, sd_2)

del_norte_standard_dev_all_61_90_winter %>% 
  kable(.,'html') %>%
  kable_styling() %>%
  scroll_box(width='250px',height='500px')
waterYear sd_2
1961 3.040037
1962 5.079023
1963 5.468558
1964 4.400891
1965 4.534426
1966 4.634512
1967 3.961144
1968 3.817035
1969 4.296435
1970 4.472217
1971 5.255746
1972 4.671240
1973 3.035805
1974 4.129045
1975 4.051736
1976 3.923459
1977 4.039750
1978 3.329846
1979 4.478687
1980 4.196422
1981 3.982774
1982 4.165896
1983 3.409845
1984 4.301257
1985 3.963049
1986 4.594538
1987 3.908011
1988 3.701619
1989 4.371704
1990 3.746830 
ggplot(del_norte_standard_dev_all_61_90_winter, aes(x = waterYear, y = sd_2))+#, color = waterYear)) +
  geom_line(size= 0.7) +
  #geom_line(aes) +
  theme_few() +
  geom_smooth(method = "lm", se=FALSE) +
  ylab('SD') + 
  xlab('Water year')

Winter season SD for water years 1961 - 1990

Mann-Kendall & Sen’s Slope

1961 - 1990 winter Mann-Kendall & Sen’s Slope

sd_mk_61_90_winter <- mk.test(del_norte_standard_dev_all_61_90_winter$sd_2)
print(sd_mk_61_90_winter)
## 
##  Mann-Kendall trend test
## 
## data:  del_norte_standard_dev_all_61_90_winter$sd_2
## z = -1.7841, n = 30, p-value = 0.07441
## alternative hypothesis: true S is not equal to 0
## sample estimates:
##            S         varS          tau 
## -101.0000000 3141.6666667   -0.2321839
sd_sens_61_90_winter <- sens.slope(del_norte_standard_dev_all_61_90_winter$sd_2)
print(sd_sens_61_90_winter)
## 
##  Sen's slope
## 
## data:  del_norte_standard_dev_all_61_90_winter$sd_2
## z = -1.7841, n = 30, p-value = 0.07441
## alternative hypothesis: true z is not equal to 0
## 95 percent confidence interval:
##  -0.043243697  0.002862457
## sample estimates:
## Sen's slope 
## -0.02253227

1991-2020

#average water year temperature

del_norte_yearly_wy_aver_91_20 <- del_norte_clean %>%
  filter(waterYear >= 1991 & waterYear <= 2020) %>% 
  group_by(waterYear) %>% 
  mutate(aver_ann_temp = mean(avg_T_c))

#Average temperature by day for all water years:

del_norte_daily_wy_aver_91_20 <- del_norte_yearly_wy_aver_91_20 %>% 
  group_by(daymonth) %>% 
  mutate(aver_day_temp = mean(avg_T_c))

#average mean temperature by day for the period of record:

del_norte_daily_wy_aver_91_20 <- del_norte_daily_wy_aver_91_20 %>% 
  group_by(daymonth) %>% 
  mutate(all_ave_temp = mean(del_norte_daily_wy_aver_91_20$aver_day_temp)) %>% 
  group_by(waterYear) %>% 
  mutate(nmbr = n())
  

del_norte_standard_dev_91_20 <- del_norte_daily_wy_aver_61_90 %>% 
  group_by(waterYear) %>% 
  mutate(residual = (all_ave_temp-aver_ann_temp)+avg_T_c-aver_day_temp) %>% 
  mutate(deviation = abs(residual-lag(residual)))

del_norte_standard_dev_91_20 <- del_norte_standard_dev_91_20 %>% 
  group_by(waterYear) %>% 
  mutate(resid_mean = mean(residual)) %>%
  mutate(sd_1 = residual-resid_mean) %>% 
  mutate(sd_2 = (((sum((sd_1)^2))/((nmbr-1))))^(0.5)) %>%
  distinct(sd_2, .keep_all = TRUE) %>% 
   select(waterYear, sd_2)


del_norte_standard_dev_91_20 %>% 
  kable(.,'html') %>%
  kable_styling() %>%
  scroll_box(width='250px',height='500px')
waterYear sd_2
1961 2.652011
1962 3.838105
1963 3.849170
1964 3.400975
1965 3.510945
1966 3.457497
1967 3.418993
1968 3.074156
1969 3.359801
1970 3.825266
1971 4.060860
1972 3.569781
1973 3.170185
1974 3.312662
1975 3.157682
1976 3.102911
1977 3.190448
1978 3.131516
1979 3.771941
1980 3.361648
1981 3.511685
1982 3.252683
1983 3.364953
1984 3.869026
1985 3.356020
1986 3.555517
1987 3.118324
1988 3.251738
1989 3.520066
1990 3.140466 
ggplot(del_norte_standard_dev_91_20, aes(x = waterYear, y = sd_2))+#, color = waterYear)) +
  geom_line(size= 0.7) +
  #geom_line(aes) +
  theme_few() +
  geom_smooth(method = "lm", se=FALSE) +
  ylab('SD') + 
  xlab('Water year')

SD for water years 1991 - 2020

1991 - 2020 Mann-Kendall & Sen’s Slope

sd_mk_91_20 <- mk.test(del_norte_standard_dev_91_20$sd_2)
print(sd_mk_91_20)
## 
##  Mann-Kendall trend test
## 
## data:  del_norte_standard_dev_91_20$sd_2
## z = -0.82069, n = 30, p-value = 0.4118
## alternative hypothesis: true S is not equal to 0
## sample estimates:
##           S        varS         tau 
##  -47.000000 3141.666667   -0.108046
sd_sens_91_20 <- sens.slope(del_norte_standard_dev_91_20$sd_2)
print(sd_sens_91_20)
## 
##  Sen's slope
## 
## data:  del_norte_standard_dev_91_20$sd_2
## z = -0.82069, n = 30, p-value = 0.4118
## alternative hypothesis: true z is not equal to 0
## 95 percent confidence interval:
##  -0.019247118  0.007185514
## sample estimates:
##  Sen's slope 
## -0.005340888

Summer and Winter 91-20

Summer temperature standard deviation

# using the 1991- 2020 data frame

del_norte_standard_dev_all_91_20_summer <- del_norte_daily_wy_aver_91_20 %>%
  filter(waterDay >= 244 & waterDay <= 335) %>% # this might be better off as daymonth rather than day of water year due to leap year
  group_by(waterYear) %>% 
  mutate(residual = (all_ave_temp-aver_ann_temp)+avg_T_c-aver_day_temp) %>% 
  mutate(deviation = abs(residual-lag(residual))) %>% 
  group_by(waterYear) %>% 
  mutate(nmbr = n())

del_norte_standard_dev_all_91_20_summer <- del_norte_standard_dev_all_91_20_summer %>% 
  group_by(waterYear) %>% 
  mutate(resid_mean = mean(residual)) %>%
  mutate(sd_1 = residual-resid_mean) %>% 
  mutate(sd_2 = (((sum((sd_1)^2))/((nmbr-1))))^(0.5)) %>%
  distinct(sd_2, .keep_all = TRUE) %>% 
   select(waterYear, sd_2)

del_norte_standard_dev_all_91_20_summer %>% 
  kable(.,'html') %>%
  kable_styling() %>%
  scroll_box(width='250px',height='500px')
waterYear sd_2
1991 1.539861
1992 1.920572
1993 1.814709
1994 1.568152
1995 2.299793
1996 1.574377
1997 1.673528
1998 1.922938
1999 1.721594
2000 1.601182
2001 1.627835
2002 1.672425
2003 2.074105
2004 2.227562
2005 1.832692
2006 1.985751
2007 1.700212
2008 1.834979
2009 1.773057
2010 1.875866
2011 1.460643
2012 1.504805
2013 2.219445
2014 1.840426
2015 1.755329
2016 1.991201
2017 1.959171
2018 1.780971
2019 1.815002
2020 1.852469 
ggplot(del_norte_standard_dev_all_91_20_summer, aes(x = waterYear, y = sd_2))+#, color = waterYear)) +
  geom_line(size= 0.7) +
  #geom_line(aes) +
  theme_few() +
  geom_smooth(method = "lm", se=FALSE) +
  ylab('SD') + 
  xlab('Water year')

Summer season SD for water years 1991 - 2020

Mann-Kendall & Sen’s Slope

Summer 91-20 standard deviations.

sd_mk_91_20_summer <- mk.test(del_norte_standard_dev_all_91_20_summer$sd_2)
print(sd_mk_91_20_summer)
## 
##  Mann-Kendall trend test
## 
## data:  del_norte_standard_dev_all_91_20_summer$sd_2
## z = 1.2489, n = 30, p-value = 0.2117
## alternative hypothesis: true S is not equal to 0
## sample estimates:
##            S         varS          tau 
##   71.0000000 3141.6666667    0.1632184
sd_sens_91_20_summer <- sens.slope(del_norte_standard_dev_all_91_20_summer$sd_2)
print(sd_sens_91_20_summer)
## 
##  Sen's slope
## 
## data:  del_norte_standard_dev_all_91_20_summer$sd_2
## z = 1.2489, n = 30, p-value = 0.2117
## alternative hypothesis: true z is not equal to 0
## 95 percent confidence interval:
##  -0.003960904  0.013198708
## sample estimates:
## Sen's slope 
## 0.006701254

Winter temperature standard deviation

# using the 91-20 data frame

del_norte_standard_dev_all_91_20_winter <- del_norte_daily_wy_aver_91_20 %>%
  filter(waterDay >= 32 & waterDay <= 182) %>% # this might be better off as daymonth rather than day of water year due to leap year
  group_by(waterYear) %>% 
  mutate(residual = (all_ave_temp-aver_ann_temp)+avg_T_c-aver_day_temp) %>% 
  mutate(deviation = abs(residual-lag(residual))) %>% 
  group_by(waterYear) %>% 
  mutate(nmbr = n())

del_norte_standard_dev_all_91_20_winter <- del_norte_standard_dev_all_91_20_winter %>% 
  group_by(waterYear) %>% 
  mutate(resid_mean = mean(residual)) %>%
  mutate(sd_1 = residual-resid_mean) %>% 
  mutate(sd_2 = (((sum((sd_1)^2))/((nmbr-1))))^(0.5)) %>%
  distinct(sd_2, .keep_all = TRUE) %>% 
   select(waterYear, sd_2)

del_norte_standard_dev_all_91_20_winter %>% 
  kable(.,'html') %>%
  kable_styling() %>%
  scroll_box(width='250px',height='500px')
waterYear sd_2
1991 4.305480
1992 3.361147
1993 3.592868
1994 3.827918
1995 4.190223
1996 4.230802
1997 4.220827
1998 3.428370
1999 3.797588
2000 3.836822
2001 3.905816
2002 3.816384
2003 3.883725
2004 4.645949
2005 4.237901
2006 4.365107
2007 4.778690
2008 4.917485
2009 3.921060
2010 3.659791
2011 4.799234
2012 3.931567
2013 5.235721
2014 4.016448
2015 4.022615
2016 3.771542
2017 3.896221
2018 3.612361
2019 3.921619
2020 3.697315 
ggplot(del_norte_standard_dev_all_91_20_winter, aes(x = waterYear, y = sd_2))+#, color = waterYear)) +
  geom_line(size= 0.7) +
  #geom_line(aes) +
  theme_few() +
  geom_smooth(method = "lm", se=FALSE) +
  ylab('SD') + 
  xlab('Water year')

Winter season SD for water years 1991 - 2020

Mann-Kendall & Sen’s Slope

1991 - 2020 winter Mann-Kendall & Sen’s Slope

sd_mk_91_20_winter <- mk.test(del_norte_standard_dev_all_91_20_winter$sd_2)
print(sd_mk_91_20_winter)
## 
##  Mann-Kendall trend test
## 
## data:  del_norte_standard_dev_all_91_20_winter$sd_2
## z = 0.67796, n = 30, p-value = 0.4978
## alternative hypothesis: true S is not equal to 0
## sample estimates:
##            S         varS          tau 
## 3.900000e+01 3.141667e+03 8.965517e-02
sd_sens_91_20_winter <- sens.slope(del_norte_standard_dev_all_91_20_winter$sd_2)
print(sd_sens_91_20_winter)
## 
##  Sen's slope
## 
## data:  del_norte_standard_dev_all_91_20_winter$sd_2
## z = 0.67796, n = 30, p-value = 0.4978
## alternative hypothesis: true z is not equal to 0
## 95 percent confidence interval:
##  -0.01252701  0.02017045
## sample estimates:
## Sen's slope 
## 0.004767841

Comparing Climate Normal means

(from 356:)

1991-2022

del_norte_yearly_wy_aver_1991_2022_CN <- del_norte_clean %>%
  filter(waterYear >= 1991 & waterYear <= 2022) %>% 
  group_by(waterYear) %>% 
  mutate(aver_ann_temp = mean(avg_T_c))

#Average temperature by day for all water years:

del_norte_daily_wy_aver_1991_2022_CN <- del_norte_yearly_wy_aver_1991_2022_CN %>% 
  group_by(daymonth) %>% 
  mutate(aver_day_temp = mean(avg_T_c))

#average mean temperature by day for the period of record:

del_norte_daily_wy_aver_1991_2022_CN <- del_norte_daily_wy_aver_1991_2022_CN %>% 
  group_by(daymonth) %>% 
  mutate(all_ave_temp = mean(del_norte_daily_wy_aver_1991_2022_CN$aver_day_temp)) %>% 
  group_by(waterYear) %>% 
  mutate(nmbr = n())
del_norte_daily_wy_aver_1991_2022_CN2 <- del_norte_daily_wy_aver_1991_2022_CN %>% 
  #filter(waterYear == "1987" | waterYear == "2021") %>%
  group_by(waterDay) %>%
  mutate(date_temp = mean(avg_T_c))%>% 
  select(waterDay, date_temp) %>% 
  distinct(waterDay, .keep_all = TRUE)
  
del_norte_daily_wy_aver_1991_2022_CN2$date_temp <- signif(del_norte_daily_wy_aver_1991_2022_CN2$date_temp,3) #reduce the sig figs

ggplot(del_norte_daily_wy_aver_1991_2022_CN2, aes(x = waterDay, y = date_temp))+#, color = waterYear)) +
  geom_line(size= 0.7) +
  theme_few()

1961-1990

del_norte_yearly_wy_aver_1961_1990_CN <- del_norte_clean %>%
  filter(waterYear >= 1961 & waterYear <= 1990) %>% 
  group_by(waterYear) %>% 
  mutate(aver_ann_temp = mean(avg_T_c))

#Average temperature by day for all water years:

del_norte_daily_wy_aver_1961_1990_CN <- del_norte_yearly_wy_aver_1961_1990_CN %>% 
  group_by(daymonth) %>% 
  mutate(aver_day_temp = mean(avg_T_c))

#average mean temperature by day for the period of record:

del_norte_daily_wy_aver_1961_1990_CN <- del_norte_daily_wy_aver_1961_1990_CN %>% 
  group_by(daymonth) %>% 
  mutate(all_ave_temp = mean(del_norte_daily_wy_aver_1961_1990_CN$aver_day_temp)) %>% 
  group_by(waterYear) %>% 
  mutate(nmbr = n())
del_norte_daily_wy_aver_1961_1990_CN2 <- del_norte_daily_wy_aver_1961_1990_CN %>% 
  #filter(waterYear == "1987" | waterYear == "2021") %>%
  group_by(waterDay) %>%
  mutate(date_temp = mean(avg_T_c))%>% 
  select(waterDay, date_temp) %>% 
  distinct(waterDay, .keep_all = TRUE)
  
del_norte_daily_wy_aver_1961_1990_CN2$date_temp <- signif(del_norte_daily_wy_aver_1961_1990_CN2$date_temp,3) #reduce the sig figs

ggplot(del_norte_daily_wy_aver_1961_1990_CN2, aes(x = waterDay, y = date_temp))+#, color = waterYear)) +
  geom_line(size= 0.7) +
  theme_few()

1931-1960

del_norte_yearly_wy_aver_1931_1960_CN <- del_norte_clean %>%
  filter(waterYear >= 1931 & waterYear <= 1960) %>% 
  group_by(waterYear) %>% 
  mutate(aver_ann_temp = mean(avg_T_c))

#Average temperature by day for all water years:

del_norte_daily_wy_aver_1931_1960_CN <- del_norte_yearly_wy_aver_1931_1960_CN %>% 
  group_by(daymonth) %>% 
  mutate(aver_day_temp = mean(avg_T_c))

#average mean temperature by day for the period of record:

del_norte_daily_wy_aver_1931_1960_CN <- del_norte_daily_wy_aver_1931_1960_CN %>% 
  group_by(daymonth) %>% 
  mutate(all_ave_temp = mean(del_norte_daily_wy_aver_1931_1960_CN$aver_day_temp)) %>% 
  group_by(waterYear) %>% 
  mutate(nmbr = n())
del_norte_daily_wy_aver_1931_1960_CN2 <- del_norte_daily_wy_aver_1931_1960_CN %>% 
  #filter(waterYear == "1987" | waterYear == "2021") %>%
  group_by(waterDay) %>%
  mutate(date_temp = mean(avg_T_c))%>% 
  select(waterDay, date_temp)  %>% 
  distinct(waterDay, .keep_all = TRUE)
  
del_norte_daily_wy_aver_1931_1960_CN2$date_temp <- signif(del_norte_daily_wy_aver_1931_1960_CN2$date_temp,3) #reduce the sig figs

ggplot(del_norte_daily_wy_aver_1931_1960_CN2, aes(x = waterDay, y = date_temp))+#, color = waterYear)) +
  geom_line(size= 0.7) +
  theme_few()

1900-1930

del_norte_yearly_wy_aver_1901_1930_CN <- del_norte_clean %>%
  filter(waterYear >= 1901 & waterYear <= 1930) %>% 
  group_by(waterYear) %>% 
  mutate(aver_ann_temp = mean(avg_T_c))

#Average temperature by day for all water years:

del_norte_daily_wy_aver_1901_1930_CN <- del_norte_yearly_wy_aver_1901_1930_CN %>% 
  group_by(daymonth) %>% 
  mutate(aver_day_temp = mean(avg_T_c))

#average mean temperature by day for the period of record:

del_norte_daily_wy_aver_1901_1930_CN <- del_norte_daily_wy_aver_1901_1930_CN %>% 
  group_by(daymonth) %>% 
  mutate(all_ave_temp = mean(del_norte_daily_wy_aver_1901_1930_CN$aver_day_temp)) %>% 
  group_by(waterYear) %>% 
  mutate(nmbr = n())
del_norte_daily_wy_aver_1901_1930_CN2 <- del_norte_daily_wy_aver_1901_1930_CN %>% 
  #filter(waterYear == "1987" | waterYear == "2021") %>%
  group_by(waterDay) %>%
  mutate(date_temp = mean(avg_T_c))%>% 
  select(waterDay, date_temp) %>% 
  distinct(waterDay, .keep_all = TRUE)


del_norte_daily_wy_aver_1901_1930_CN2$date_temp <- signif(del_norte_daily_wy_aver_1901_1930_CN2$date_temp,3) #reduce the sig figs

ggplot(del_norte_daily_wy_aver_1901_1930_CN2, aes(x = waterDay, y = date_temp))+#, color = waterYear)) +
  geom_line(size= 0.7) +
  theme_few()

climate_norm <- left_join(del_norte_daily_wy_aver_1991_2022_CN2, del_norte_daily_wy_aver_1961_1990_CN2, by= 'waterDay') 

climate_norm <- left_join(climate_norm, del_norte_daily_wy_aver_1931_1960_CN2, by= 'waterDay')

climate_norm <- left_join(climate_norm, del_norte_daily_wy_aver_1901_1930_CN2, by= 'waterDay')

ggplot(climate_norm, aes(x=waterDay)) +  
         #scale_size_manual(values=c("FALSE"=5,"TRUE"=8)) +
  #scale_color_manual(values=c("FALSE"='#CC0000',"TRUE"='black')) +
geom_line(aes(y = date_temp.x), size =1, color = "orange")+
  geom_line(aes(y = date_temp.y), size =1, color = "darkgreen")+
  geom_line(aes(y = date_temp.x.x), size =1, color = "purple")+
  geom_line(aes(y = date_temp.y.y), size =1, color= "red")+
  theme_few()+
  scale_colour_manual("", breaks=c("date_temp.x","date_temp.y","date_temp.x.x", "date_temp.y.y"), values=c("date_temp.x" = "orange", "date_temp.y" = "darkgreen", "date_temp.x.x" = "purple", "date_temp.y.y" = "red"))+
  xlab("Day of Year")+
  ylab("Average temperature")